External validation, recalibration, and clinical utility of the prognostic model kidney failure risk equation in patients with CKD stages G3-4: a nationwide retrospective cohort analysis in Peru

Main Analysis - Winsorization 1.5%: 4. External Validation

Author

Percy Soto Becerra

1 Setup

rm(list = ls())

# Use pacman to check whether packages are installed, if not load
if (!require("pacman")) install.packages("pacman")
library(pacman)

# Unload all package to begin in a session with only base packages
pacman::p_unload(all)

# Install packages
pacman::p_load(
  rio, 
  here, 
  tidyverse, 
  gtools, 
  kableExtra, 
  patchwork, 
  flextable, 
  furrr, 
  parallel, 
  mice, 
  survival,
  prodlim, 
  cmprsk,
  riskRegression,
  splines
)

source(here("Code", "source", "kfre_pi.R"))
source(here("Code", "source", "kfre_pr.R"))
source(here("Code", "source", "oe_ratio.R"))
source(here("Code", "source", "calibration_intercept.R"))
source(here("Code", "source", "calibration_slope.R"))
source(here("Code", "source", "auc.R"))
source(here("Code", "source", "auc_brier_boot.R"))
source(here("Code", "source", "validate.mids.R"))
source(here("Code", "source", "pool.validate.mids.R"))
source(here("Code", "source", "pool.auc.mids.R"))
source(here("Code", "source", "process_imp_cal_plot.R"))
source(here("Code", "source", "predict.mira.R"))
source(here("Code", "source", "performance_measures.R"))
source(here("Code", "source", "tidy_performance_stack.R"))    
source(here("Code", "source", "tidy_pool.R")) 
source(here("Code", "source", "process_imp_cal_plot2.R")) 

seed <- 2014

2 Load data

# # Import data
imp.datosA <- readRDS(here::here("Data", "Tidy", "Main-Winsorize-1_5", "data_impA.rds")) 
gc()
            used  (Mb) gc trigger   (Mb)  max used  (Mb)
Ncells   2832971 151.3    5252691  280.6   4101325 219.1
Vcells 117037097 893.0  151470922 1155.7 117048872 893.1

3 Configurate cores for parallel processing

n_cores <- detectCores()
n_cores <- (n_cores - 2) / 2
# Evaluate futures in parallel - max of two workers to avoid hogging resources
future::plan("multisession", workers = n_cores)

4 Pre-processing data

imp.datosA <- imp.datosA |> 
  mutate(risk2y = kfre_pr(imp.datosA, 2),
         risk5y = kfre_pr(imp.datosA, 5), 
         eventdf = factor(eventd)) |> 
  select(.imp, .id, time, eventd, eventdf, risk2y, risk5y) |>  
  filter(.imp != 0)

5 Predictive Performance

5.1 Calibration and Discrimination Measures

future::plan("multisession", workers = n_cores)
results_stack3a4_2y <- tidy_performance_stack(imp.datosA, 
                                  horizon = 2, 
                                  primary_event = 1, 
                                  pred = "risk2y",
                                  seed = seed, 
                                  n_cores = n_cores)

rio::export(results_stack3a4_2y , here("Data", "Tidy", "Main-Winsorize-1_5", 
                                       "res_valext__kfre_orig_stack3a4_2y.rds"))

future::plan("multisession", workers = n_cores)
results_stack3a4_5y <- tidy_performance_stack(imp.datosA, 
                                  horizon = 5, 
                                  primary_event = 1, 
                                  pred = "risk5y",
                                  seed = seed, 
                                  n_cores = n_cores)

rio::export(results_stack3a4_5y, here("Data", "Tidy", "Main-Winsorize-1_5", 
                                      "res_valext_kfre_orig_stack3a4_5y.rds"))
res_pool1 <- tidy_pool(results_stack3a4_2y) 

res_pool1 |> 
  kbl() |> 
  kable_classic(full_width = T, html_font = "Cambria")
term estimate ll ul pval ubar b t dfcom df riv lambda fmi n
Average predicted risk 0.0148486 NA NA NA NaN 0.0000001 NaN 30030 NA NaN NaN NaN 30031
Overall observerd risk 0.0273046 0.0254293 0.0291799 0.0000000 0.0000009 0.0000000 0.0000009 30030 30030.0000 0.0000000 0.0000000 0.0000666 30031
Log OE ratio 0.6092915 0.5329617 0.6856213 0.0000000 0.0012279 0.0002845 0.0015152 30030 2473.0857 0.2340161 0.1896378 0.0007311 30031
OE difference 0.0124560 0.0105161 0.0143959 0.0000000 0.0000009 0.0000001 0.0000010 30030 12686.3357 0.0699669 0.0653917 0.0001525 30031
Calibration Intercept 0.0248474 -0.1324316 0.1821265 0.7565925 0.0043575 0.0020444 0.0064224 30030 914.7292 0.4738640 0.3215113 0.0018290 30031
Calibration Slope 0.5824941 0.5329646 0.6320236 0.0000000 0.0005467 0.0000906 0.0006382 30029 4054.5128 0.1674280 0.1434161 0.0004576 30031
Logit AUC 1.9605067 1.8436399 2.0773735 0.0000000 0.0030191 0.0005287 0.0035531 30030 3740.4228 0.1768742 0.1502915 0.0004941 30031
Brier 0.0247911 0.0231400 0.0264421 0.0000000 0.0000007 0.0000000 0.0000007 30030 25245.4258 0.0239194 0.0233606 0.0000783 30031
res_pool2 <- tidy_pool(results_stack3a4_5y) 

res_pool2 |> 
  kbl() |> 
  kable_classic(full_width = T, html_font = "Cambria")
term estimate ll ul pval ubar b t dfcom df riv lambda fmi n
Average predicted risk 0.0447519 NA NA NA NaN 0.0000005 NaN 30030 NA NaN NaN NaN 30031
Overall observerd risk 0.0476187 0.0450853 0.0501521 0.0000000 0.0000017 0.0000000 0.0000017 30030 30030.0000 0.0000000 0.0000000 0.0000666 30031
Log OE ratio 0.0622124 0.0007751 0.1236496 0.0471823 0.0007368 0.0002418 0.0009810 30030 1491.6060 0.3314888 0.2489610 0.0011716 30031
OE difference 0.0028668 -0.0000193 0.0057529 0.0515475 0.0000017 0.0000005 0.0000022 30030 1753.1433 0.2961194 0.2284662 0.0010088 30031
Calibration Intercept -0.4678735 -0.5889528 -0.3467942 0.0000000 0.0023383 0.0014492 0.0038020 30030 644.6382 0.6259778 0.3849854 0.0024937 30031
Calibration Slope 0.5801802 0.5423074 0.6180529 0.0000000 0.0003212 0.0000515 0.0003732 30029 4260.6666 0.1618482 0.1393024 0.0004364 30031
Logit AUC 1.7741859 1.6870515 1.8613202 0.0000000 0.0018599 0.0001150 0.0019761 30030 14229.1899 0.0624473 0.0587769 0.0001364 30031
Brier 0.0419170 0.0397886 0.0440455 0.0000000 0.0000011 0.0000001 0.0000012 30030 11823.2705 0.0746902 0.0694993 0.0001632 30031
tab_res_2y <- res_pool1 |> 
  select(term, estimate, ll, ul) |> 
  mutate(
    across(c(estimate, ll, ul), ~ if_else(term == "Log OE ratio", exp(.x), .x)), 
    across(c(estimate, ll, ul), ~ if_else(term == "Logit AUC", exp(.x) / (1  + exp(.x)), .x)), 
    term = if_else(term == "Log OE ratio", "OE ratio", term), 
    term = if_else(term == "Logit AUC", "AUC", term), 
    across(c(estimate, ll, ul), ~ if_else(term %in% 
                                            c("Average predicted risk", 
                                              "Overall observerd risk", 
                                              "OE difference"), 100 * .x, .x)), 
    across(c(estimate, ll, ul), ~ round(.x, 2)), 
    measures = case_when(term == "Average predicted risk" ~ str_glue("{estimate}%"), 
                       term %in% c("Overall observerd risk", "OE difference") ~ str_glue("{estimate}% ({ll}% to {ul}%)"),
                       TRUE ~ str_glue("{estimate} ({ll} to {ul})")
                       )
    ) |> 
  select(term, measures) |> 
  mutate(grupo = c(rep("Calibration", 6), "Discrimination", "Overall performance"), 
         year = "2 years")

tab_res_5y <- res_pool2 |> 
  select(term, estimate, ll, ul) |> 
  mutate(
    across(c(estimate, ll, ul), ~ if_else(term == "Log OE ratio", exp(.x), .x)), 
    across(c(estimate, ll, ul), ~ if_else(term == "Logit AUC", exp(.x) / (1  + exp(.x)), .x)), 
    term = if_else(term == "Log OE ratio", "OE ratio", term), 
    term = if_else(term == "Logit AUC", "AUC", term), 
    across(c(estimate, ll, ul), ~ if_else(term %in% 
                                            c("Average predicted risk", 
                                              "Overall observerd risk", 
                                              "OE difference"), 100 * .x, .x)), 
    across(c(estimate, ll, ul), ~ round(.x, 2)), 
    measures = case_when(term == "Average predicted risk" ~ str_glue("{estimate}%"), 
                       term %in% c("Overall observerd risk", "OE difference") ~ str_glue("{estimate}% ({ll}% to {ul}%)"),
                       TRUE ~ str_glue("{estimate} ({ll} to {ul})")
                       )
    ) |> 
  select(term, measures) |> 
  mutate(grupo = c(rep("Calibration", 6), "Discrimination", "Overall performance"), 
         year = "5 years")

tab_res0 <- tab_res_2y |>
  bind_rows(tab_res_5y)

tab_res0 |> 
  as_grouped_data(groups = "year") |> 
  as_grouped_data(groups = "grupo") |> 
  flextable::as_flextable(hide_grouplabel = TRUE) |> 
  autofit() |> 
  set_header_labels(
    year = "Time horizon", 
    term = "Performance measure", 
    measures = "Original KFRE"
  ) |>  
  bold(j = 1) |> 
  set_caption("Table. Performance measures of KFRE in the external dataset of patients with CKD Stages 3a-4") |>  
  theme_booktabs() |>   
  bold(bold = TRUE, part = "header") -> table_perf_final

table_perf_final %>% 
  flextable::save_as_docx(path = here("Tables", "Main-Winsorize-1_5", "Table_Imputed_Performance_Modelo_Original.docx"))

table_perf_final

Time horizon

Performance measure

Original KFRE

2 years

Calibration

Average predicted risk

1.48%

Overall observerd risk

2.73% (2.54% to 2.92%)

OE ratio

1.84 (1.7 to 1.99)

OE difference

1.25% (1.05% to 1.44%)

Calibration Intercept

0.02 (-0.13 to 0.18)

Calibration Slope

0.58 (0.53 to 0.63)

Discrimination

AUC

0.88 (0.86 to 0.89)

Overall performance

Brier

0.02 (0.02 to 0.03)

5 years

Calibration

Average predicted risk

4.48%

Overall observerd risk

4.76% (4.51% to 5.02%)

OE ratio

1.06 (1 to 1.13)

OE difference

0.29% (0% to 0.58%)

Calibration Intercept

-0.47 (-0.59 to -0.35)

Calibration Slope

0.58 (0.54 to 0.62)

Discrimination

AUC

0.85 (0.84 to 0.87)

Overall performance

Brier

0.04 (0.04 to 0.04)

rm(list=ls()[! ls() %in% c("imp.datosA","vdata", "primary_event", "horizon", 
                           "process_imp_cal_plot", "seed", "n_cores")])
gc()
           used  (Mb) gc trigger  (Mb)  max used   (Mb)
Ncells  3177787 169.8    5252691 280.6   5252691  280.6
Vcells 24374928 186.0   93104695 710.4 181150158 1382.1

5.2 Moderate calibration: Calibration curves lowess based on subdistributional hazards

primary_event <- 1

n_internal_knots <- 5

# Seleccion del grupo: Stages 3-4----

# A 2 años----
horizon <- 2

vdata <- imp.datosA %>% 
  rename(pred = risk2y) |> 
  select(.imp, .id, time, eventd, pred) |> 
  mutate(cll_pred = log(-log(1 - pred)))

rcs_vdata <- ns(vdata$cll_pred, df = n_internal_knots + 1)
colnames(rcs_vdata) <- paste0("basisf_", colnames(rcs_vdata))
knots <- attr(rcs_vdata, "knots")
bound.knots <-  attr(rcs_vdata, "Boundary.knots")

pp <- seq(min(vdata$pred), max(vdata$pred), length.out = 1000)
cll_pp <- log(-log(1 - pp))
rcs_cll_pp <- ns(cll_pp, knots = knots, Boundary.knots = bound.knots)
colnames(rcs_cll_pp) <- paste0("basisf_", colnames(rcs_cll_pp))

vdata_bis_pp <- cbind(pp, as.data.frame(rcs_cll_pp))

future::plan("multisession", workers = n_cores)
subdist_df_imp <- future_map(1:max(vdata$.imp),
                         process_imp_cal_plot, 
                         vdata = vdata, 
                         primary_event = primary_event, 
                         horizon = horizon, 
                         type = "subdist_hazard", 
                         n_internal_knots = n_internal_knots, 
                         vdata_bis_pp, 
                         .options = furrr_options(seed = seed, 
                                                  packages = c("riskRegression", 
                                                               "survival", 
                                                               "splines", 
                                                               "cmprsk",
                                                               "tidyverse")), 
                         .progress = TRUE)
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
# 5 knots seems to give somewhat equivalent graph to pseudo method with bw = 0.05
subdist_df_imp_obs <- subdist_df_imp |> 
  bind_rows() |> 
  filter(type == "observed")

subdist_df_stack <- subdist_df_imp_obs |>
  group_by(.imp) |> 
  mutate(deciles_risk = as.integer(quantcut(risk, seq(0, 1, by = 0.1)))) |> 
  group_by(.imp, deciles_risk) |> 
  summarise(obs_mean_imp = mean(obs), 
            risk_mean_imp = mean(risk)) |> 
  group_by(deciles_risk) |> 
  summarise(obs_mean_pool = mean(obs_mean_imp), 
            risk_mean_pool = mean(risk_mean_imp))
  
subdist_df_imp_fix <- subdist_df_imp |> 
  bind_rows() |>  
  filter(type == "fixed") |> 
  arrange(risk) |> 
  summarise(obs_pool = mean(obs), 
            .by = risk) |> 
  mutate(deciles_risk = quantcut(risk, seq(0, 1, by = 0.1)))

rio::export(subdist_df_imp, here("Data", "Tidy", "Main-Winsorize-1_5", "subdist_df_imput_3a4_2y.rds"))
rio::export(subdist_df_stack, here("Data", "Tidy", "Main-Winsorize-1_5", "subdist_df_deciles_3a4_2y.rds"))

cal_mod2y_imp <- subdist_df_imp_obs |> 
  select(.imp, risk, obs) |> 
  summarise(ICI = mean(abs(risk - obs)), 
            E50 = quantile(abs(risk - obs), c(0.5)), 
            E90 = quantile(abs(risk - obs), c(0.9)), 
            Emax = max(abs(risk - obs)), 
            root_square_bias = sqrt(mean((risk - obs) ^ 2)), 
            .by = ".imp"
  ) 

cal_mod2y_pool <- cal_mod2y_imp |> 
  summarise(ICI = mean(ICI), 
            E50 = mean(E50), 
            E90 = mean(E90), 
            Emax = mean(Emax), 
            root_square_bias = mean(root_square_bias)) |> 
  mutate(year = "2 years") |> 
  select(year, everything())

rio::export(cal_mod2y_imp, here("Data", "Tidy", "Main-Winsorize-1_5", "cal_Original2y_imp.rds"))

# Grafico de calibracion
ggplot() +
  geom_abline(intercept = 0, slope = 1, colour = "red", linetype = 2) + 
  geom_line(data = subdist_df_imp_fix, 
            aes(x = risk, y = obs_pool),
            alpha = 0.5, color = "black") +
  geom_point(data = subdist_df_stack,
             aes(x = risk_mean_pool, y = obs_mean_pool),
             shape = 23,
             stroke = 0.1,
             fill = "blue") + 
  scale_y_continuous(breaks = seq(0, 1, 0.2), limits = c(0, 1)) +
  scale_x_continuous(breaks = seq(0, 1, 0.2), limits = c(0, 1)) +
  theme_bw() + 
  labs(x = "Predicted risks", y = "Observed outcome proportions") + 
  # geom_rug(data = subdist_df_stack,
  #          aes(x = risk_mean_pool, y = obs_mean_pool),
  #          alpha = 0.1) +
  coord_fixed(ratio = 1, expand = TRUE)  -> plot_mod0_2y

# Grafico de calibracion con curvas imputadas
plot_mod0_2y + 
  geom_line(data = subdist_df_imp_obs, 
            aes(x = risk, y = obs, group = .imp),
            alpha = 0.1, color = "#38B8F7") -> plot_mod0_imp_2y

# A 5 años----
horizon <- 5

vdata <- imp.datosA %>% 
  rename(pred = risk5y) |> 
  select(.imp, .id, time, eventd, pred) |> 
  mutate(cll_pred = log(-log(1 - pred))) 

rcs_vdata <- ns(vdata$cll_pred, df = n_internal_knots + 1)
colnames(rcs_vdata) <- paste0("basisf_", colnames(rcs_vdata))
knots <- attr(rcs_vdata, "knots")
bound.knots <-  attr(rcs_vdata, "Boundary.knots")

pp <- seq(min(vdata$pred), max(vdata$pred), length.out = 1000)
cll_pp <- log(-log(1 - pp))
rcs_cll_pp <- ns(cll_pp, knots = knots, Boundary.knots = bound.knots)
colnames(rcs_cll_pp) <- paste0("basisf_", colnames(rcs_cll_pp))

vdata_bis_pp <- cbind(pp, as.data.frame(rcs_cll_pp))

future::plan("multisession", workers = n_cores)
subdist_df_imp <- future_map(1:max(vdata$.imp),
                         process_imp_cal_plot, 
                         vdata = vdata, 
                         primary_event = primary_event, 
                         horizon = horizon, 
                         type = "subdist_hazard", 
                         n_internal_knots = n_internal_knots, 
                         vdata_bis_pp, 
                         .options = furrr_options(seed = seed, 
                                                  packages = c("riskRegression", 
                                                               "survival", 
                                                               "splines", 
                                                               "cmprsk",
                                                               "tidyverse")), 
                         .progress = TRUE)
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
[1] "basisf_1" "basisf_2" "basisf_3" "basisf_4" "basisf_5" "basisf_6"
# 5 knots seems to give somewhat equivalent graph to pseudo method with bw = 0.05
subdist_df_imp_obs <- subdist_df_imp |> 
  bind_rows() |> 
  filter(type == "observed")

subdist_df_stack <- subdist_df_imp_obs |>
  group_by(.imp) |> 
  mutate(deciles_risk = as.integer(quantcut(risk, seq(0, 1, by = 0.1)))) |> 
  group_by(.imp, deciles_risk) |> 
  summarise(obs_mean_imp = mean(obs), 
            risk_mean_imp = mean(risk)) |> 
  group_by(deciles_risk) |> 
  summarise(obs_mean_pool = mean(obs_mean_imp), 
            risk_mean_pool = mean(risk_mean_imp))
  
subdist_df_imp_fix <- subdist_df_imp |> 
  bind_rows() |>  
  filter(type == "fixed") |> 
  arrange(risk) |> 
  summarise(obs_pool = mean(obs), 
            .by = risk) |> 
  mutate(deciles_risk = quantcut(risk, seq(0, 1, by = 0.1)))

rio::export(subdist_df_imp, here("Data", "Tidy", "Main-Winsorize-1_5", "subdist_df_imput_3a4_5y.rds"))
rio::export(subdist_df_stack, here("Data", "Tidy", "Main-Winsorize-1_5", "subdist_df_deciles_3a4_5y.rds"))

cal_mod5y_imp <- subdist_df_imp_obs |> 
  select(.imp, risk, obs) |> 
  summarise(ICI = mean(abs(risk - obs)), 
            E50 = quantile(abs(risk - obs), c(0.5)), 
            E90 = quantile(abs(risk - obs), c(0.9)), 
            Emax = max(abs(risk - obs)), 
            root_square_bias = sqrt(mean((risk - obs) ^ 2)), 
            .by = ".imp"
  ) 

cal_mod5y_pool <- cal_mod5y_imp |> 
  summarise(ICI = mean(ICI), 
            E50 = mean(E50), 
            E90 = mean(E90), 
            Emax = mean(Emax), 
            root_square_bias = mean(root_square_bias)) |> 
  mutate(year = "5 years") |> 
  select(year, everything())

rio::export(cal_mod5y_imp, here("Data", "Tidy", "Main-Winsorize-1_5", "cal_Original5y_imp.rds"))
# Grafico de calibracion
ggplot() +
  geom_abline(intercept = 0, slope = 1, colour = "red", linetype = 2) + 
  geom_line(data = subdist_df_imp_fix, 
            aes(x = risk, y = obs_pool),
            alpha = 0.5, color = "black") +
  geom_point(data = subdist_df_stack,
             aes(x = risk_mean_pool, y = obs_mean_pool),
             shape = 23,
             stroke = 0.1,
             fill = "blue", 
             alpha = 0.5) + 
  scale_y_continuous(breaks = seq(0, 1, 0.2), limits = c(0, 1)) +
  scale_x_continuous(breaks = seq(0, 1, 0.2), limits = c(0, 1)) +
  theme_bw() + 
  labs(x = "Predicted risks", y = "Observed outcome proportions") + 
  # geom_rug(data = subdist_df_stack,
  #          aes(x = risk_mean_pool, y = obs_mean_pool),
  #          alpha = 0.1) +
  coord_fixed(ratio = 1, expand = TRUE)  -> plot_mod0_5y

# Grafico de calibracion con curvas imputadas
plot_mod0_5y + 
  geom_line(data = subdist_df_imp_obs, 
            aes(x = risk, y = obs, group = .imp),
            alpha = 0.3, color = "#38B8F7") -> plot_mod0_imp_5y

## Grafico final
plot_mod0_2y <- plot_mod0_2y + 
  labs(title = "Original Model\n(2 year KFRE)") + 
  theme(plot.title = element_text(hjust = 0.5))

plot_mod0_5y <- plot_mod0_5y + 
  labs(title = "Original Model\n(5 year KFRE)") + 
  theme(plot.title = element_text(hjust = 0.5))

plot_mod0_imp_2y <- plot_mod0_imp_2y +
  labs(title = "Original Model\n(2 year KFRE)") +
  theme(plot.title = element_text(hjust = 0.5))

plot_mod0_imp_5y <- plot_mod0_imp_5y +
  labs(title = "Original Model\n(5 year KFRE)") +
  theme(plot.title = element_text(hjust = 0.5))

plot_cal_mod0 <- (plot_mod0_2y | plot_mod0_5y) + plot_annotation(tag_levels = "a")
plot_cal_imp_mod0 <- (plot_mod0_imp_2y | plot_mod0_imp_5y) + plot_annotation(tag_levels = "a")

export(plot_mod0_2y, here("Data", "Tidy", "Main-Winsorize-1_5", "plot_mod0_2y.rds"))
export(plot_mod0_5y, here("Data", "Tidy", "Main-Winsorize-1_5", "plot_mod0_5y.rds"))
export(plot_cal_mod0, here("Data", "Tidy", "Main-Winsorize-1_5", "plot_cal_mod0.rds"))

export(plot_mod0_imp_2y, here("Data", "Tidy", "Main-Winsorize-1_5", "plot_mod0_imp_2y.rds"))
export(plot_mod0_imp_5y, here("Data", "Tidy", "Main-Winsorize-1_5", "plot_mod0_imp_5y.rds"))
export(plot_cal_imp_mod0, here("Data", "Tidy", "Main-Winsorize-1_5", "plot_cal_imp_mod0.rds"))


ggsave(filename = "Plot_Calibration_imputed.jpeg", 
       device = "jpeg", 
       plot = plot_cal_mod0, 
       path = here("Figures", "Main-Winsorize-1_5"), 
       scale = 2, 
       width = 12, 
       height = 6,
       units = "cm", 
       dpi = 1000)

ggsave(filename = "Plot_Calibration_imputed_estabilidad.jpeg", 
       device = "jpeg", 
       plot = plot_cal_imp_mod0, 
       path = here("Figures", "Main-Winsorize-1_5"), 
       scale = 2, 
       width = 12, 
       height = 6,
       units = "cm", 
       dpi = 1000)

gc()
           used  (Mb) gc trigger   (Mb)  max used   (Mb)
Ncells  3343159 178.6    5252691  280.6   5252691  280.6
Vcells 85954162 655.8  272778960 2081.2 272766285 2081.1
knitr::include_graphics(here("Figures", "Main-Winsorize-1_5", "Plot_Calibration_imputed.jpeg"))

Calibration curves for each group and prediction horizon
knitr::include_graphics(here("Figures", "Main-Winsorize-1_5", "Plot_Calibration_imputed_estabilidad.jpeg"))

Calibration curves for each group and prediction horizon

5.3 Moderate Calibration

cal_mod2y_pool |> 
  bind_rows(cal_mod5y_pool) |> 
  pivot_longer(cols = ICI:root_square_bias, 
               names_to = "metrica", 
               values_to = "modelo") |> 
  mutate(modelo = round(modelo, 3)) |> 
  as_grouped_data(groups = "year") |> 
  flextable::as_flextable(hide_grouplabel = TRUE) |> 
  autofit() |>  
  set_header_labels(
    year = "Time horizon", 
    metrica = "Performance measure", 
    modelo = "Original KFRE"
  ) |> 
  bold(i = c(1, 7), j = 1) |> 
  set_caption("Table. Moderate Calibration Metrics") |>  
  theme_booktabs() |>   
  bold(bold = TRUE, part = "header") -> table_cal_mod

table_cal_mod %>% 
  flextable::save_as_docx(path = here("Tables", "Main-Winsorize-1_5", "Table_Imputed_Performance_Moderate_Cal_Modelo_Original.docx"))

table_cal_mod

Performance measure

Original KFRE

2 years

ICI

0.016

E50

0.005

E90

0.049

Emax

0.541

root_square_bias

0.031

5 years

ICI

0.017

E50

0.005

E90

0.035

Emax

0.475

root_square_bias

0.049

6 Reproducibility Ticket

sessionInfo()
R version 4.4.1 (2024-06-14)
Platform: x86_64-pc-linux-gnu
Running under: Ubuntu 24.04 LTS

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.12.0 
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.12.0

Random number generation:
 RNG:     L'Ecuyer-CMRG 
 Normal:  Inversion 
 Sample:  Rejection 
 
locale:
 [1] LC_CTYPE=es_ES.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=es_ES.UTF-8        LC_COLLATE=es_ES.UTF-8    
 [5] LC_MONETARY=es_ES.UTF-8    LC_MESSAGES=es_ES.UTF-8   
 [7] LC_PAPER=es_ES.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=es_ES.UTF-8 LC_IDENTIFICATION=C       

time zone: America/Lima
tzcode source: system (glibc)

attached base packages:
[1] splines   parallel  stats     graphics  grDevices utils     datasets 
[8] methods   base     

other attached packages:
 [1] riskRegression_2023.12.21 cmprsk_2.2-12            
 [3] prodlim_2024.06.25        survival_3.7-0           
 [5] mice_3.16.0               furrr_0.3.1              
 [7] future_1.34.0             flextable_0.9.6          
 [9] patchwork_1.2.0           kableExtra_1.4.0         
[11] gtools_3.9.5              lubridate_1.9.3          
[13] forcats_1.0.0             stringr_1.5.1            
[15] dplyr_1.1.4               purrr_1.0.2              
[17] readr_2.1.5               tidyr_1.3.1              
[19] tibble_3.2.1              ggplot2_3.5.1            
[21] tidyverse_2.0.0           here_1.0.1               
[23] rio_1.2.2                

loaded via a namespace (and not attached):
  [1] rstudioapi_0.16.0       jsonlite_1.8.8          shape_1.4.6.1          
  [4] magrittr_2.0.3          TH.data_1.1-2           jomo_2.7-6             
  [7] farver_2.1.2            nloptr_2.1.1            rmarkdown_2.28         
 [10] geepack_1.3.11          ragg_1.3.2              vctrs_0.6.5            
 [13] minqa_1.2.8             askpass_1.2.0           base64enc_0.1-3        
 [16] htmltools_0.5.8.1       polspline_1.1.25        broom_1.0.6            
 [19] Formula_1.2-5           mitml_0.4-5             parallelly_1.38.0      
 [22] htmlwidgets_1.6.4       sandwich_3.1-0          zoo_1.8-12             
 [25] uuid_1.2-1              lifecycle_1.0.4         iterators_1.0.14       
 [28] pkgconfig_2.0.3         Matrix_1.7-0            R6_2.5.1               
 [31] fastmap_1.2.0           digest_0.6.37           numDeriv_2016.8-1.1    
 [34] colorspace_2.1-1        rprojroot_2.0.4         textshaping_0.4.0      
 [37] Hmisc_5.1-3             fansi_1.0.6             timechange_0.3.0       
 [40] compiler_4.4.1          fontquiver_0.2.1        withr_3.0.1            
 [43] htmlTable_2.4.3         backports_1.5.0         highr_0.11             
 [46] R.utils_2.12.3          pan_1.9                 MASS_7.3-61            
 [49] lava_1.8.0              quantreg_5.98           openssl_2.2.1          
 [52] tools_4.4.1             foreign_0.8-86          zip_2.3.1              
 [55] future.apply_1.11.2     nnet_7.3-19             R.oo_1.26.0            
 [58] glue_1.7.0              mets_1.3.4              nlme_3.1-165           
 [61] grid_4.4.1              checkmate_2.3.2         cluster_2.1.6          
 [64] generics_0.1.3          gtable_0.3.5            tzdb_0.4.0             
 [67] R.methodsS3_1.8.2       data.table_1.16.0       hms_1.1.3              
 [70] xml2_1.3.6              utf8_1.2.4              foreach_1.5.2          
 [73] pillar_1.9.0            lattice_0.22-5          SparseM_1.84-2         
 [76] tidyselect_1.2.1        fontLiberation_0.1.0    rms_6.8-2              
 [79] knitr_1.48              fontBitstreamVera_0.1.1 gridExtra_2.3          
 [82] svglite_2.1.3           xfun_0.47               stringi_1.8.4          
 [85] yaml_2.3.10             pacman_0.5.1            boot_1.3-30            
 [88] evaluate_0.24.0         codetools_0.2-20        officer_0.6.6          
 [91] gdtools_0.4.0           cli_3.6.3               rpart_4.1.23           
 [94] systemfonts_1.1.0       munsell_0.5.1           Rcpp_1.0.13            
 [97] globals_0.16.3          MatrixModels_0.5-3      jpeg_0.1-10            
[100] lme4_1.1-35.5           listenv_0.9.1           glmnet_4.1-8           
[103] viridisLite_0.4.2       mvtnorm_1.3-1           timereg_2.0.6          
[106] scales_1.3.0            rlang_1.1.4             multcomp_1.4-26